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Single-particle excitations under coexisting electron correlation and disorder: a 
numerical study of the Anderson-Hubbard model 
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Interplay of electron correlation and randomness is studied by using the Anderson-Hubbard 
model within the Hartree-Fock approximation. Under the coexistence of short-range interaction 
and diagonal disorder, we obtain the ground-state phase diagram in three dimensions, which 
includes an antiferromagnetic insulator, an antiferromagnetic metal, a paramagnetic insulator 
(Anderson-localized insulator) and a paramagnetic metal. Although only the short-range in- 
teraction is present in this model, we find unconventional soft gaps in the insulating phases 
irrespective of electron filling, spatial dimensions and long-range order, where the single-particle 
density of states (DOS) vanishes with a power-law scaling in one dimension (ID) or even faster 
in two dimensions (2D) and three dimensions (3D) toward the Fermi energy. We call it soft 
Hubbard gap. Moreover, exact-diagonalization results in ID support the formation of the soft 
Hubbard gap beyond the mean-field level. The formation of the soft Hubbard gap cannot be 
attributed to a conventional theory by Efros and Shklovskii (ES) owing the emergence of soft 
gaps to the long-range Coulomb interaction. Indeed, based on a picture of multivalley energy 
landscape, we propose a phenomenological scaling theory, which predicts a scaling of the DOS, 
A in energy E as A(E) oc exp[— (— 7 log \E — Ep\) ]. Here, d is the spatial dimension, Ef is the 
Fermi energy and 7 is a non-universal constant. This scaling is in perfect agreement with the 
numerical results. We further discuss a correction of the scaling of the DOS by the long-range 
part of the Coulomb interaction, which modifies the scaling of Efros and Shklovskii. Further- 
more, explicit formulae for the temperature dependence of the DC resistivity via variable-range 
hopping under the influence of the soft gaps are derived. Finally, we compare the present theory 
with experimental results of SrRui-^Ti^Os. 

KEYWORDS: electron correlation, disorder, Anderson-Hubbard model, single-particle density of states, 
soft gap, variable-range hopping 



1. Introduction 

1.1 Single-particle gaps in insulators 

Strongly-correlated electron systems continue to be a 
challenging issue of the condensed matter physics. Es- 
pecially, metal-insulator transitions have been attracting 
much attention, because of various phases found in their 
vicinities. 1 When the electron correlation becomes dom- 
inant compared to the kinetic energy, the ground state 
undergoes a transition from a metal into a correlation- 
induced insulator. The Mott transition at specific elec- 
tron density is a typical example. 2 In the Mott insula- 
tor, the electrons are localized on the individual atomic 
orbitals to avoid the strong on-site Coulomb repulsion, 
leading to the opening of a gap in the single-particle DOS 
as illustrated in Fig. 1 (a). There only the spin and or- 
bital degrees of freedom remain at low energies. The ab- 
sence of the single-particle excitations at low energies is 
a common feature of correlation-induced insulators, such 
as antiferromagnetic insulators and charge-ordered insu- 
lators. 

Another source of the localization in the strongly- 
correlated electron systems is disorder (or randomness), 
which is inevitably present in real materials. The disor- 
der drives the metal-insulator transition as the Anderson 
transition. 3-5 There, the insulators arc characterized not 
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(a) Insula tor in duced by inte ractio n 

(b) Anderson Insulator 




Fig. 1. (Color online) Schematic of insulators with different types 
of gaps in the single-particle density of states: (a) Insulator in- 
duced by interaction, (b) Anderson insulator, (c) Disorder depen- 
dence of a single-particle gap and (d) Insulator with soft gap. 



by the vanishing carrier number but by a vanishing relax- 
ation time accompanying the quantum localization of the 
wave functions by the impurity scattering. In contrast to 
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the Mott insulator, the Anderson insulators without elec- 
tron correlations exhibit no gap as illustrated in Fig. 1 
(b), indicating that the gapless single-particle excitations 
are essential in determining their physical properties at 
low energies. This makes the Anderson insulators com- 
pletely different from the Mott insulators. 

Since the electron correlation and randomness in- 
evitably coexist in real materials, clarification of the 
single-particle excitations under the coexistence is impor- 
tant for understanding their physical properties. How- 
ever, since the Mott and Anderson insulators have qual- 
itatively different low-energy excitations, the DOS at 
low energies under their coexistence is highly-non-trivial, 
which is the main topic of this paper. Let us consider how 
a single-particle gap behaves when fluctuating random 
potentials are introduced as disorder into the pure Mott 
insulator. When the disorder strength is weak enough, 
the single-particle gap might survive despite appearance 
of localized impurity levels induced within the gap as 
illustrated in Fig. 1 (c). With the increase of disorder 
strength W, width of distribution of random potentials, 
however, the gap amplitude gradually decreases and fi- 
nally the gap collapses at W = Wqap- One might think 
that for W > Wgap, the gap completely closes and the 
DOS at Ef becomes nonzero. Namely, the ground state 
undergoes a transition to a simple Anderson insulator, 
where the gapless single-particle excitations dominate 
the low-energy physics. Unfortunately, this naive expec- 
tation is not correct. 

In their seminal work, Efros and Shklovskii 6 (ES) con- 
sidered an amorphous or a doped crystalline semiconduc- 
tor, where the Coulomb interaction is nearly unscreened. 
They showed that assuming nonzero A(Ep) in the ground 
state in the presence of the long-range Coulomb inter- 
action, the ground state is unstable against a particle- 
hole excitation between localized states at Ef, which 
is not a single-particle excitation but a multiple exci- 
tation from the ground state. In other words, this means 
that the supposed ground state with nonzero A(Ef) is 
relaxed by electron-hole excitations to reduce A(Ef) un- 
der the influence of the long-range Coulomb interaction. 
As a result of this excitonic effect, soft Coulomb gap in- 
deed opens and the DOS at Ep remains zero even for 
W > Wgap as illustrated in Fig. 1 (d), where the DOS 
is scaled 

A(E) oc IS-SfI^ 1 . (1.1) 

Here d is the spatial dimension. Note that although we 
focus on half filling in Fig. 1, the soft Coulomb gap is 
not restricted to half filling. 

It is important to emphasize that this power-law scal- 
ing is not restricted to the critical point which sepa- 
rates the correlation-induced and Anderson insulators. 
Instead, according to the ES theory, this power-law scal- 
ing generically dominates the whole insulating phase, 
which indicates that the insulator region is always criti- 
cal. Indeed, the existence of the Coulomb gap was con- 
firmed by numerical simulations 8-10 and in electron tun- 
neling experiments 11-13 later. 



1.2 Recent studies on coexistence of short-range inter- 
action and disorder 

When the screening gets stronger e.g., near metal- 
insulator transitions with divergence of the dielectric 
constant, the effects of the long-range part of the 
Coulomb interaction are restricted to lower and lower 
energies. In this case, the soft Coulomb gap arising from 
the ES mechanism is expected to shrink to extremely 
low energies and the effect of the short-range part prac- 
tically determines electronic structures in the experimen- 
tally accessible energy scale. However, within the ES the- 
ory, short-range interactions do not generate soft gaps, 
because the excitonic effect is negligible for the short- 
range interaction. Therefore, according to the ES theory, 
one could speculate that the soft gap vanishes in the in- 
sulator near the metal-insulator transition. 

Although the coexistence of the interaction and ran- 
domness is common in strongly-correlated materials, par- 
ticularly perovskite-type compounds with B-site sub- 
stitution ABi_2,B^,03 offer a suitable stage for the in- 
vestigation of the combined effects. SrRui_j;Ti x 03 is a 
promising candidate for this purpose. One of the end 
component SrRuOs is a correlated ferromagnetic metal 
(Tq = 165 K); the other end component SrTiOs is a 
band insulator with a wide band gap (~ 3.2 eV). The Ad 
Ru band is located at the Fermi energy, while the Ti 3d 
band is well separated from the Ru Ad band around the 
Fermi energy. Thus the Ti atoms act as impurities for 
the itinerant Ru 4c? electrons, and the potential height 
of the introduced disorder is as much as of the order 
of 1 eV. As a consequence, around 0.3 < x < 0.5, the 
material undergoes a metal-insulator transition into a 
correlated Anderson insulator. 14 In the vicinity of the 
metal-insulator transition, indeed, recent photocmission 
results of SrRui_ x Tij;03 15 ' 16 indicates breakdown of the 
ES scaling in 3D. Although the power law with the ex- 
ponent a = 2 is expected from the ES theory for 3D, 
the exponent obtained by the fitting of the experimental 
data (~ 1.2) is clearly different from a = 2. The break- 
down of the ES theory has been observed in other mate- 
rials such as LaNii-^Mn^Os, 17 indicating the ubiquity 
of unconventional soft gaps in the vicinities of the metal- 
insulator transitions. The deviation from the ES theory 
near the metal-insulator transition was also observed in a 
numerical study with the unrestricted Hartree-Fock ap- 
proximation. 18 

Despite long-time efforts to go beyond the ES theory, 
however, effects of short-range interaction have not been 
fully understood yet because of the difficulties in han- 
dling quantum effects, which is important for short-range 
interaction. The Anderson-Hubbard model is one of the 
minimal models of real strongly-correlated materials un- 
der the coexistence of the interaction and randomness. 
The Anderson-Hubbard Hamiltonian is defined by 

with N s sites and N e electrons, where t is a hopping in- 
tegral, U is the on-site repulsion, c\ a (c ia ) is the creation 
(annihilation) operator for an electron with spin a on the 
site i. The number operator is defined by = c\ a Ci a 
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and \i is the chemical potential. In addition to the usual 
Hubbard Hamiltonian representing the itinerancy and 
the short-range interaction (namely, the first and the sec- 
ond terms of eq. (1.2)), disorder is represented by the 
spatially uncorrelated spin-independent random poten- 
tial Vi. 

Many numerical techniques have been applied to the 
Anderson-Hubbard model. For example, the quantum 
Mote Carlo (QMC) method was applied in ID, 19 ' 20 2D 21 
and 3D. 22 Furthermore, the dynamical mean-field the- 
ory was extended to disordered systems. 23 At infinite di- 
mensions, Dobrosavljcvic and Kotliar formulated a vari- 
ant of the dynamical mean-field theory (DMFT), 23,24 so- 
called statistical DMFT, which is exact in the limit of the 
infinite coordination number or in the non-interacting 
limit. By using the site-dependent bath functions, sta- 
tistical DMFT can partially treat spatial correlations. 
Subsequently, Dobrosavljevic et al. derived a mean-field 
theory, 25 which is simpler but ignores the spatial cor- 
relations. Another approach is the (unrestricted) site- 
dependent Hartrce-Fock (HF) approximation. Although 
it treats the electron correlation in the mean-field level, 
it is certainly beyond the mean-field theory of the type 
that allows only spatially uniform mean fields. In fact, it 
can describe the inhomogeneity of the electronic struc- 
tures by using site-dependent mean- fields. Tusch et al. 
obtained the ground-state phase diagram including both 
of the magnetic and charge degrees of freedom in 3D. 26 

Recently, several numerical studies have reported sup- 
pression of the DOS at the Fermi energy even for short- 
range interaction. For example, soft gaps were reported 
in a HF study in 3D. 27 Although they claimed a power- 
law scaling of the DOS with a ~ 0.5, the origin of the soft 
gap has not been clarified at all. Recent numerical studies 
in 2D by using the exact diagonalization showed a dip of 
the DOS near Ep. 28 These strongly suggest the presence 
of an unconventional mechanism which suppresses the 
DOS even with short-range interaction. However, a nu- 
merical study with statistical DMFT 23 claimed nonzero 
A(Ep) in the insulating phase, and even the divergence 
of A(E-p) toward the metal- insulator transition from the 
metallic side, which completely disagrees with the HF 
results. On the other hand, in several mean-field studies 
which ignores the spatial correlations, no singularity was 
reported in the DOS. 25,30 We clearly need further stud- 
ies for comprehensive understanding of the short-range 
case. 

In this paper, through numerical analyses of the 3D 
Anderson-Hubbard model, we show that there exists 
a soft gap even though only short-range interaction is 
present in the Anderson-Hubbard model. We call this 
unconventional soft gap soft Hubbard gap. We show nu- 
merical evidences of the soft Hubbard gaps within the 
HF approximation in ID, 2D and 3D. Further support 
by the exact diagonalization in ID is given. In order to 
clarify the origin of the soft Hubbard gap, we propose a 
phenomenological theory based on a picture of multival- 
ley energy landscape, which corresponds to emergence of 
many excited states degenerated with the ground state. 
Our scaling theory predicts an unconventional scaling 
of the density of states A(E) in energy E as A(E) cx 



exp[— (— 7 log \E — Ep\) d }. Here 7 is a non-universal con- 
stant. We show that this predicted scaling is consistent 
with the numerically-observed scaling. A part of the nu- 
merical evidences for the soft Hubbard gap and the scal- 
ing theory have already been given in a letter briefly. 31 In 
this paper, however, we analyze the ground-state phase 
diagram of the 3D Anderson-Hubbard model in greater 
detail. Especially, the criticality of the metal-insulator 
transitions is clarified from a viewpoint of the formation 
of the soft Hubbard gap. Furthermore, by extending our 
scaling theory, we clarify effects of the long-range part of 
the Coulomb interaction responsible for low-energy exci- 
tations. Especially, we show that the ES theory is seri- 
ously modified in the presence of the long-range Coulomb 
interaction, when we consider the multiply-excited states 
by extending our scaling theory originally constructed 
for the short-range interaction. In order to inspire fu- 
ture experimental efforts to examine the validity of the 
present fundamental proposal, we study temperature de- 
pendence of the DC resistivity in the presence of the 
soft gap. We further compare the experimental results 
for SrRui_a,Ti x 03 with the present theory. 

This paper is organized as follows; In § 2, we intro- 
duce the Anderson-Hubbard model and numerical meth- 
ods employed. In § 3, we show numerical results of 
the Anderson-Hubbard model in 3D and ID within the 
Hartree-Fock approximation as well as by the exact di- 
agonalization. Section 4 is devoted to the scaling theory 
of the soft gap and its extension to discrete distribu- 
tion functions of random potentials and the long-range 
Coulomb interaction. In § 5, we derive transport prop- 
erties in the presence of the soft gap. Comparisons with 
experimental results are also given. The summary and 
discussion are given in § 6. 

2. Model and Method 

2.1 Anderson-Hubbard Model 

In this paper, we analyze the Anderson-Hubbard 
model, whose Hamiltonian is defined by eq. (1.2). We 
employ a cubic lattice for d = 3, a square lattice for 
d = 2 and a chain lattice for d = 1. We take the lattice 
spacing as the length unit. The spin-independent random 
potential Vi representing randomness is assumed to fol- 
low two models of the distribution Pv(Vi): the box type 
of width 2W, 



\Vi\ < W) 
(otherwise) 



1 

2ir 



(2.1) 



with the average (Vi) = 0, and the Gaussian type, 

1 ( Yl 

2a 2 



iV(Vi) = 



2ira 



e*P "7^2 > ( 2 -2) 



where a 2 = W 2 /12 with the average (Vi) — 0. Because 
W is proportional to the standard deviation of the dis- 
tribution for the both models, W is the parameter to 
control the strength of disorder. For both the two distri- 
butions, (i = U/2 corresponds to half filling. 

2.2 Hartree-Fock approximation 

We first employ the Hartree-Fock (HF) approxima- 
tion, where the wave function is approximated by a sin- 
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gle Slater determinant consisting of a set of orthonor- 
mal single-particle orbitals {4> n } (n is an orbital index). 
Within this trial wave function, the variational principle 
leads to the HF equation as 

{Ho + U y^X( n il) n i1 + ( n it) n il)}<t>n = £n4>n, (2-3) 
i 

where TIq is the one-body part of the Hamiltonian and we 
neglect (c^Cii). Here (ni a ) arc the site-dependent mean 
fields. Later, we will show that the inhomogeneity of the 
electronic structures is necessary for the formation of the 
soft Hubbard gap. To find a site-dependent mean-field 
solution (ni a ) for the HF equations, we employ the iter- 
ative scheme starting from an appropriate initial guess 
until the convergence condition, |An i(T | < 10~ 5 (Vic) is 
satisfied. Here Ani a denotes a change in the mean field 
at the site i with spin ct, rn a before and after an iteration. 
Initial guesses of mean fields employed in our calculations 
will be described later in each case. In order to accelerate 
the convergence of the iteration, we employ the Ander- 
son mixing, 32 which is easily applicable to the 3D case 
because its computational cost scales as 0(N S ) and it is 
not heavy. 

In general, many stable mean-field solutions may co- 
exist for a given realization of the random potentials. 
Therefore, after repeating the calculations for several dif- 
ferent initial guesses of the mean fields as the ground 
state, we employ the solution that has the lowest energy. 

2.3 Exact diagonalization 

In this paper, we further employ the exact diagonal- 
ization method in ID. Because the exact diagonalization 
method takes into account quantum fluctuations ignored 
by the Hartree-Fock approximation, it is suitable for the 
examination of robustness of our Hartree-Fock results 
against the quantum fluctuations. By diagonalizing the 
full Hamiltonian matrix by using LAPACK routines, 33 
we obtain all the eigenstates. Then the single-particle 
density of states is calculated by 

£„E 4 J(^c + l|ct|0;iV c )| 2 

A{E) = { xKE-E + {n)) (E > E F ) 

E„EJW-l|4|0;iV e >| 

x5(E-E-(n)) (E < E F ), 

where N e is the number of electrons in the ground state, 
\n; N) is the n-th eigenstates with N electrons (n > 0), 
and |0; N c ) is the ground state. The single-particle exci- 
tation energies E±(n) arc defined as 

E±(n) =E- E n (N c ± 1) + E (N C ) - E F . (2.5) 

Here E n (N) is the energy of the n-th eigenstates with N 
electrons, |0; N c ). In order to reduce the computational 
cost, we divide the full Hilbcrt space into subspaces with 
a fixed set of (iVj , JVj.), where JVf and N± denote the num- 
ber of up-spin electrons and that of down-spin electrons, 
respectively. The maximum dimension of the subspaces 
is 400 (N s = 6). 

3. Numerical Results 

In this section, we analyze the Anderson-Hubbard 
model in ID and 3D numerically. 



Our main new results are the following: 

(1) Determination of the ground-state phase diagram of 
the Anderson-Hubbard model in 3D within the HF 
approximation, 

(2) Discovery of the unconventional soft gap driven by 
the short-range interaction regardless of electron fill- 
ing and spatial dimensionality, 

(3) Clarification of the unconventional scaling of the soft 
gap; the DOS decays faster than any power law in 
3D, while it follows a power law in ID, 

(4) Clarification of the criticality of the metal-insulator 
transitions. 

We first analyze the Anderson-Hubbard model in 3D 
within the HF approximation. Further numerical evi- 
dences of the soft gap in ID are given within the Hartree- 
Fock approximation. We also show cxact-diagonalization 
results in ID beyond the mean-field level. Detailed scal- 
ing analyses within the Hartree-Fock approximation in 
2D will be given in the next section § 4.1. In the fol- 
lowing, we focus on half filling unless otherwise stated. 
Although a part of the results have already been reported 
briefly, 31 in this paper, we show more detailed numerical 
analyses which were omitted in the previous letter, e.g. 
the criticality of the metal- insulator transitions. 

3.1 Three dimensions 

3.1.1 Ground- state phase diagram 

We study the ground state of the 3D Anderson- 
Hubbard model with the Gaussian distribution of Py 
in detail. A part of the results have already been re- 
ported briefly. 31 In our 3D study, we take the hopping 
integral t as the energy unit. Figure 2 shows the cal- 
culated ground-state phase diagram within the HF ap- 
proximation. At U = 0, the Anderson-Hubbard model 
undergoes a metal-insulator transition (Anderson transi- 
tion) from the paramagnetic metal (PM) to the param- 
agnetic insulator (PI) at a finite strength of the disorder, 
W c = 21.29 ± 0.02. 34 On the other hand, at W = 0, 
since the system is half-filled with a perfect nesting, the 
ground state is the antiferromagnetic insulator (AFI) for 
any nonzero value of U. Here, we discuss the ground- 
state phase diagram for U, W > 0. First, we focus on the 
spin degree of freedom. For W > 0, the ground state is 
paramagnetic near U = 0. With increasing the interac- 
tion, the ground states undergoes an antiferromagnetic 
transition at a critical point U c (> 0). Within the reso- 
lution of our calculation, U c monotonically increases as 
the disorder strength W increases. Next, we focus on the 
charge degree of freedom. The ground state is insulating 
for U, W 1, which contains AFI as well as paramag- 
netic insulator (PI) (PI is usually identified as Anderson 
insulator). Metallic phases are restricted to a dome- like 
region (U < 6 and W < 25). In addition to a paramag- 
netic metal (PM), we found an antiferromagnetic metal 
(AFM). In the following, we will show detailed analyses 
of the ground-state phase diagram. 
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W c = 21.29 ±0.02 

A 



5 10 15 20 25 30 W 



Fig. 2. (Color online) Ground-state phase diagram of 3D 
Anderson-Hubbard model at half filling for Gaussian distribu- 
tion of Py. Abbreviations are: AFI, antiferromagnetic insula- 
tor; AFM, antiferromagnetic metal; PI, paramagnetic insula- 
tor (Anderson insulator); PM, paramagnetic metal. The dot- 
ted line denotes the asymptotic scaling of the antiferromagnetic 
transition for U S> t and W t. The critical value of the 
metal-insulator transition in the non-interacting limit is given 
by W c = 21.29 ± 0.02. 34 



3.1.2 Antiferromagnetic transitions 

Wc first discuss the antiferromagnetic transition. The 
antiferromagnetic order parameter m(Q) is defined by 



m(Q) = 



B 5 * 



JQr, 



(3.1) 



where Q = (ir, ir, ir). At W = 0, since the system is half- 
filled with the perfect nesting condition satisfied, the an- 
tiferromagnetic order vanishes in an essentially-singular 
way toward U = as 39 



t ( 2wt 
m (Q) = jj exp I — — 



(3.2) 



While the amplitude of the antiferromagnetic gap given 
by Um(Q) also follows essentially-singular scaling and is 
small near U = 0, U c increases sharply from zero with 
the increasing disorder W. 

We consider the antiferromagnetic transition in the 
atomic limit (U, W ^> t). In this limit, there are two types 
of sites: singly occupied sites for \Vi\ < U/2, and doubly 
occupied or empty sites for V% < —U/2 or Vi > U/2, re- 
spectively. Thus the Anderson-Hubbard model reduces 
to a site-diluted antiferromagnet, 40 where the antifer- 
romagnetic transition coincides with percolation of the 
singly occupied sites within the mean-field level, because 
of the antiferromagnetic interaction between neighboring 
singly occupied sites. Thus U c is determined by 
»+t/ c 

P(Vi)dVi = Pc , (3.3) 

where p c ~ 0.311608 is the site-percolation threshold of 
the cubic lattice. By solving this equation for the Gaus- 
sian distribution eq. (2.2), one obtains U c /W ~ 0.2315 
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Fig. 3. (Color online) Square of antiferromagnetic order parame- 
ter m(Q) 2 at U = 3, W = 20 and W = 30 for (a), (b) and (c), 
respectively. The solid lines are the fit by the mean-field critical- 
ly p = 1/2. Results for 4 X 4 X 4 and 6x6x6 are shown only 
for the open-shell boundary condition. In (b) and (c), we only 
show results from the open-shell boundary condition. 



in the atomic limit. This asymptotic scaling is plotted as 
the dotted line in Fig. 2. As we will show below, the phase 
boundary numerically obtained starts following this scal- 
ing with increasing W and U . 

In order to determine the antiferromagnetic transition 
line in the intermediate region, we calculate the anti- 
ferromagnetic order parameter with cubic unit cells of 
L x L x L. We employ two boundary conditions: closed- 
shell and open-shell boundary conditions. At U = and 
W = 0, the periodic boundary condition corresponds 
to open-shell configurations for L = An and closed-shell 
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configurations for L = An + 2, where n is a positive inte- 
ger. In contrast, the antiperiodic boundary condition in 
3D corresponds to closed-shell configurations for L = An 
and open-shell configurations for L = 4n + 2. Therefore 
we define the open-shell boundary condition as the pe- 
riodic boundary condition for L = An and the antiperi- 
odic boundary condition for L = 4n + 2. The closed-shell 
boundary condition is defined as the opposite cases. For 
finite-size lattices, the antiferromagnetic order parame- 
ter calculated with the open-shell boundary condition is 
larger than that with the closed-shell boundary condi- 
tion, because of an enhancement of the electron correla- 
tion for the open-shell boundary condition. 

We employ four kinds of initial guesses of mean fields: 

(1) Uniform and antifcrromagnctically-ordered states, 

(2) Ground states in the limit of t — ► degenerating 
with each other within the second-order perturba- 
tion with respect to t, which consists of locally- 
antifcrromagnctically-ordercd clusters parted by 
doublons and holons, 

(3) Paramagnetic states with uniform charge distribu- 
tion, 

(4) States with random charge-spin distributions. 

All these initial guesses are further perturbed by dif- 
ferent choices of random noise in both of the spin and 
charge sectors before putting into the self-consistent it- 
eration. We typically need several tens of initial guesses 
for convergence of the antiferromagnetic order param- 
eters. Figure 3 (a) shows the calculated antiferromag- 
netic order parameter at U = 3 with the open-shell and 
closed-shell boundary conditions, respectively. Error bars 
smaller than the symbols are dropped for the simplicity 
of the figure. For each boundary condition, as the sys- 
tem size increases, the antiferromagnetic order param- 
eter well converges to the mean-field critical behavior 
of j3 = 1/2. We estimate the upper and lower limits of 
the critical point by fitting of the data with the open- 
shell boundary condition and those with the closed-shell 
boundary condition, respectively. 

At W/t = 20, 30, we show the antiferromagnetic order 
parameter calculated only with the open-shell bound- 
ary condition, because difference in the calculated an- 
tiferromagnetic order parameter was found to be negli- 
gible between the two boundary conditions. As shown 
in Fig. 3 (b) and (c), the antiferromagnetic order pa- 
rameter again well converges to the mean-field critical 
behavior, as the system size becomes larger. The esti- 
mated antiferromagnetic transition points of U c /W = 
0.231 ± 0.003, 0.217 ± 0.003 for W/t = 20,30 are close 
to the atomic-limit value (U c /W ~ 0.2315), establishing 
the validity of our calculations. The error bars are de- 
termined by the fitting of the data with the system size 
12 x 12 x 12 for W = 20, 10 x 10 x 10 for W = 30, 
respectively. This correct asymptotic behavior of U c was 
not reproduced in the previous HF study. 26 

On the other hand, K. Byczuk et al. found a non- 
monotonic behavior of U c for the antiferromagnetic tran- 
sition as a function of the increasing disorder strength W 
within DMFT with spin degrees of freedom. 35 Namely, 
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Fig. 4. (Color online) Inverse of spin structure factor S(Q) mul- 
tiplied by the number of sites N B as a function of U at W = 20. 



they found that U c becomes zero as a function of W 
at the Anderson transition point in the non-interacting 
limit. In contrast, we found a monotonic increase of U c 
in our 3D study, which we believe more reasonable. We 
show a further numerical evidence at W = 20. We define 
spin structure factor S(q) as 



S(q) 



1 



Sj)e 



iq-iri-Tj) 



*,3 



(3.4) 



In Fig. 4, we plot spin structure factor at Q — (tt, 7T, 7t) 
for W = 20. The quantity N S S(Q) diverges toward the 
antiferromagnetic transition point approaching from the 
limit of U — 0, which further supports the existence of 
the antiferromagnetic transition around U/t ~ 4. 

3.1.3 Metal-insulator transitions 

Next, we discuss the metal-insulator transitions. We 
identify insulating phases by extrapolation of the local- 
ization lengths £ to the bulk limit. The localization length 
£ is defined by the asymptotic behavior of single-particle 
orbitals near E-p at long distances as 



cx exp(-r/£), 



(3.5) 



where r is the distance from the center of the orbital. 
In order to observe the asymptotic behavior easily, we 
employ pseudo-ID unit cells of Lx Lx M, where M L 
as has been employed in Refs. 36-38. In fact, we adopt 
M = 1000, 250, 250, 100 for L = 4,6, 8, 10, respectively. 
We use single-particle orbitals within an energy window 
of the width 0.01 around the Fermi energy, namely, \E — 
Ep\ < 0.01 to calculate the localization lengths. 

Figure 5 shows the extrapolation of the calculated lo- 
calization lengths for U = 2,4,5,6 to the bulk limit. In 
insulating phases, localization lengths are extrapolated 
to finite values. For U < 6 and with increasing W from 
0, metals appear from AFI as in the 2D result, 40 with 
further reentrant transition to insulators (AFI or PI) at 
larger W. Even in antifcrromagnctically-ordered phases, 
there exist metal-like regions, where the antiferromag- 
netic order parameter is locally reduced due to relatively 
large fluctuations of Vi. When W is small enough, the 
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system remains insulating even though the DOS is gap- 
less, because these metal-like regions are small in size and 
isolated from each other. With increasing W, however, 
they grow in size and finally percolates, corresponding 
to the first metal-insulator transition into metals (AFM 
or PM). Further increase of the disorder i.e., W ^> t, 
results in the reentranst metal-insulator transition into 
insulators (AFI or PI). It should be noted that because 
the Gaussian-type distribution is unbounded, the gaplcss 
regions appear and the antifcrromagnctic gap closes as 
soon as the disorder strength W becomes nonzero. This 
indicates Wqap = 0. Thus there is no real gap for W > 
in the case of the Gaussian- type distribution of Py. How- 
ever, this is not a generic feature for Py bounded in a 
finite region. 

In the previous DMFT study, K. Byczuk et al. found 
a first-order transition between the Mott insulator and 
the metal. 30 In contrast, the double occupancy in our 
3D study exhibits a jump neither at the metal-insulator 
transitions nor at the antiferromagnetic transitions as 
shown in Fig. 6, indicating that they are not of the first- 
order but continuous. Now, we discuss the criticality of 
the metal- insulator transitions. In the non-interacting 
limit, the normalized localization lengths Al = £l/L for 
the different system sizes cross at the Anderson transi- 
tion. Namely, A^ is constant with respect to the width 
of the ID bar L when L is large enough. Furthermore, 
the critical parameter A c = Al(W c ) was reported to 
be universal, not depending on the distribution of the 
random potentials. 34 In Fig. 7, we show scaling plots of 
the localization lengths at U = 2, U = 4 and U = 5. 
Even for the interacting cases, the normalized localiza- 
tion lengths indeed seem to cross at the metal-insulator 
transitions with a universal critical parameter A c ~ 0.8. 
This value is clearly different from that of the non- 
interacting limit (A c = 0.576 ± 0.002 34 ), indicating that 
the metal-insulator transitions for U > belong to a uni- 
versality class different from that of the Anderson tran- 
sition at U = 0. Further support for this unconventional 
criticality will be given in § 3.1.4 in the context of the 
density of states. 

Before moving on the next section, we analyze the en- 
ergy dependence of the localization length in the para- 
magnetic insulator (PI). In Fig. 8 (a), we show the lo- 
calization length extrapolated to the bulk limit, £ at 
t = 1, U = 4 and W = 30 as well as at t = 1, 
C/ = 4 and W = 30. Although the localization length 
is nearly independent of the energy near the Fermi en- 
ergy at U = 0, it depends on the energy in a more 
complex manner at U = 4. Around the Fermi energy, 
namely \E — Ep\ < 0.5, the inverse localization length 
l/£ is extrapolated to a finite value, being consistent 
with the fact that the ground state is insulating. Also 
for \E — E-p\ > 5, the extrapolated localization length is 
finite and increases toward high energies. In the inter- 
mediate region, namely 0.5 < \E — E-p\ < 5, however, 
l/£ is extrapolated to zero, indicating the existence of 
extended states. This reentrant transition in the imme- 
diate vicinity of the Fermi energy is clearly beyond the 
conventional picture of the Anderson transitions without 
electron correlations, where there is no reentrant transi- 
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Fig. 5. (Color online) Quasi-ID localization lengths as a function 
of inverse length l/L. 



tion with respect to the energy. 47 

In order to reinforce the results obtained by the 
analyses of the localization lengths, we further analyze 
another order parameter for the Anderson transition, 
the geometrically-averaged DOS, The geometrically- 
averaged DOS, Aq(E) is defined as 



A G (E) = exp[— J2^g(Ai(E)) 



where the local DOS at the site i is given by 



A(E) = --1mG iti (E + iri), 



(3.6) 



(3.7) 



which is exact in the limit of 77 = +0. Here 77 (> 
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Fig. 6. (Color online) Double occupancy with the open-shell 
boundary condition as functions of U. The system sizes are 
8x8x8 for W = 5 and 20, 10 X 10 X 10 for the rest. The 
double occupancy increases as W increases for a fixed value of 
U. This is consistent with the fact that the double occupancy 
has its maximum value 0.5 in the limit W = +co, where only an 
equal number of doubly occupied sites and empty sites exist in 
the ground state. It should be mentioned that the double occu- 
pancy is 0.25 at U = and W = 0, and in the limit U = +oo. 



0) is a broadening parameter. In the bulk limit, the 
geometrically-averaged DOS at the Fermi energy is an- 
other order parameter for the Anderson transitions, be- 
ing nonzero when states at the Fermi energy are extended 
and zero when they become localized. This is because 
when wave functions at a given energy become local- 
ized in real space, the local DOS at that energy becomes 
discrete: localized wave functions at the same energy arc 
well separated in position so the local DOS at that energy 
goes to zero far from the wave function centers. In Fig. 8 
(b), we show the geometrically-averaged DOS extrapo- 
lated to the bulk limit and further to the limit of r\ = +0. 
At low but nonzero energies, the geometrically-averaged 
DOS is extrapolated to a nonzero value, indicating the 
existence of extended states. This is consistent with the 
results obtained by the extrapolation of the localization 
length. We will discuss the implication of this unusual 
feature in the next subsection. 

3.1.4 Soft Hubbard gaps 

Next, we discuss the DOS. Figure 9(a) shows the 
DOS for typical parameters. Naively one might expect 
A(Ep) > for W > 0. Indeed, there exists no soft gap in 
the DOS at the parameter A depicted in Fig. 2 which cor- 
responds to the non-interacting Anderson insulator. The 
DOS A(E-p) is also nonzero for t = as is seen in the case 
D in Fig. 2. However, we find a soft gap over the entire in- 
sulating phases in the case of U > regardless of the an- 
tiferromagnetic order. Figure 9(a) shows the DOS for the 
typical parameters B and C, respectively. These parame- 
ters are chosen to be located in the insulating phases far 
from the metal-insulator transition, because the energy 
scale of the soft gap decreases toward the metal-insulator 
transition, which makes numerical analyses of the DOS 
difficult. A similar soft gap is seen in lower dimensions. 
Actually, one- and two-dimensional data will be analyzed 
in § 3.2 and § 4.1, respectively. We call this unconven- 
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Fig. 7. (Color online) Finite-size scaling plot of localization 
lengths at U = 2, U = 4 and (7 = 5. The metal-insulator transi- 
tion points are denoted by the centers of the circles. 



tional soft gap soft Hubbard gap, since it is driven by the 
short-range interaction U. The significance of the soft 
Hubbard gap is clear because the soft gap is established 
irrespective of the spatial dimension and electron filling. 
It should be mentioned that the soft Hubbard gap is not 
restricted to half filling as well, as we will see in ID. Al- 
though the soft gap is observed generically, its formation 
certainly has to satisfy minimal requirement: Not only 
the coexistence of the interaction and randomness but 
also the itinerancy is required for their formation. In- 
deed, the soft Hubbard gap vanishes at the parameter D 
by switching off the transfer. All of the three terms in the 
Anderson-Hubbard Hamiltonian (eq. (1.2)), namely itin- 
erancy, interaction and randomness are imperative for 
the formation of the soft gap. 
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Geometrically-averaged DOS extrapolated to the limit of r\ = +0 
and the bulk limit. The shade denotes the energy region where 
the states are extended (Aq = 0) 



This soft gap moreover bears generically an uncon- 
ventional nature: Although a power law scaling A(E) oc 
\E — E-p\ a with exponents 0.5 < a < 1 looks fit in the 
range \E — E-p\ > 0.1 being consistent with the previous 
HF study, 27 a closer look for \E — E-p\ < 0.1 shows clear 
deviation from the power-law scaling, and even faster 
decay of the DOS near Ep as shown in Fig. 9 (b). The 
deviation from a simple power law strongly suggests that 
the soft Hubbard gap originates from a novel mechanism. 

The formation of the soft Hubbard gap in the insulat- 
ing side further supports an unconventional universality 
class of the metal-insulator transitions for the interact- 
ing case U > (see § 3.1.3). The critical exponent of the 
density of states /3 is defined by 

A(E F ) = \W -W c f (Metallic side). (3.8) 

For U — 0, the critical exponent (3 is zero because the 
density of states remains nonzero in the Anderson insula- 
tor. However, the critical exponent (3 should be nonzero 
for U > because of the formation of the soft Hub- 
bard gap. This further supports that the metal-insulator 
transition for U > belongs to a universality class dif- 
ferent from that for U = 0. The formation of the soft 
Hubbard gap in the insulating phases may be respon- 
sible for the increase of the critical parameter A c from 
the non-interacting value (A c = 0.576 ± 0.002) to ~ 0.8 
by switching on the electron correlation. Now we con- 
sider how localized states at the Fermi energy become 



extended toward the metal-insulator transition by ap- 
proaching from the strongly-localized limit. In the pres- 
ence of the soft Hubbard gap, the mean distance between 
localized states near the Fermi energy for U > is larger 
than that at U — 0. This prevents the localized states 
from being hybridized with each other and percolating. 
Thus the single-particle wave functions at the Fermi en- 
ergy remain localized in the bulk limit for U > 0, even if 
the quantity A^ = L has reached the non-interacting 
critical parameter A c . 

Furthermore, the formation of the soft Hubbard is re- 
sponsible for the low-lying extended excited states in the 
insulating phases (see Fig. 8). Because of the formation 
of the soft Hubbard gap at the Fermi energy, and the 
weight excluded and transferred from the low-energy re- 
gion around the Fermi energy, the DOS forms a peak 
right outside the soft gap as seen in the DOS at U = 4 
and W = 30 (Fig. 9). Since the mean distance between lo- 
calized states in an energy window of the width AE cen- 
tered at the energy E is proportional to (A(E)AE)~ 1 / d , 
the mean distance between localized states is shortest at 
the peak position. Thus, with decreasing W, the local- 
ized states at the peak position become extended before 
those at other energies, leading to the appearance of the 
low-lying extended excited states in the vicinity of the 
metal-insulator transition. Therefore the reentrant local- 
ization and the low-lying extended excited states may be 
ubiquitous in insulators with a soft gap. 

3.2 One dimension 

In order to clarify how the unconventional scaling ob- 
served in 3D depends on the spatial dimensionality and 
is modified in lower dimensions, we further numerically 
analyze the Anderson- Hubbard model in ID with the box 
distribution of Py in detail here. 

3.2.1 Hartree-Fock approximation 

We find the unconventional soft gaps also in ID re- 
gardless of electron filling within the HF approximation. 
Figure 10 (a) shows the DOS for the hole-doped case 
as well as for the half-filled case. We employ the peri- 
odic boundary condition. Here, holes are doped with the 
chemical potential [i shifted by —10/3 from the half fill- 
ing to increase the average distance between electrons to 
capture long-range asymptotic behavior easily. Surpris- 
ingly in contrast to the 3D case, they fit well with a power 
law A{E) oc \E — E-p\ a even at low energies as shown in 
Fig. 10 (b) regardless of electron filling. Although the 
Efros-Shklovskii theory predicts A(Ep) > in ID (see 
cq. (1.1)), a logarithmic law of the soft gap was derived 
in ID by considering excitations with multi electron-hole 
pairs as 41 



A(E) 



1 



(3.9) 



log(E /\E - E \Y 

However, this conventional theory cannot explain the ob- 
served power law, which clearly indicates the existence of 
the unconventional mechanism of the soft gap. Although 
the gaps again shrink with the decreasing energy scale 
when t or U becomes smaller, the power law still holds 
when t and U are nonzero as shown in Fig. 11. 
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Fig. 9. (Color online) DOS with system size 8 x 8 x 250: A (t = 
l t U = 0,W = 30), B (t = 1,U = 4,W = 30), C (t = 1,17 = 
6, = 5), D (t = 0,U = 4, W = 30). (a) Linear plot, (b) 
Double logarithmic plot. The solid lines denote the fitting lines 
by a power law. 

3.2.2 Exact diagonalization 

In order to investigate effects of quantum fluctuations 
not taken into account in the HF approximation, we fur- 
ther analyze the Anderson-Hubbard model by using the 
exact diagonalization in ID. In Fig. 12, we further show 
the calculated DOS at t = 1, U = 10 and W = 10. 
Because we need all the eigenstates to calculate the ex- 
citation spectra of each ensemble with a fixed configu- 
ration of randomness and the ensemble average over an 
extensive number of random configurations is required, 
we restrict the system size N s to N s < 6. We find a 
dip in the DOS at Ep, where the DOS decreases as the 
system size N s increases. This indicates the formation 
of a soft Hubbard gap beyond the mean-field level. In 
order to clarify the scaling of the DOS, we employ a 
finite-size scaling analysis. We assume a scaling func- 
tion, A(e, A^T 1 ) = N-^f(eN^ /a ) = e a g(N^ /a e- 1 ) cor- 
responding to A(e, /V" 1 = 0)(Xf" and A(e = 0, TV" 1 ) oc 
N- fj (e=\E- E F \). As shown in Fig. 12, the DOS well 
converges to this scaling function with a = 0.075 and 
[3 = 0.375. Although a possible logarithmic scaling can- 
not be excluded because of the limitation in the system 
sizes of the present exact-diagonalization calculation, the 
exact-diagonalization results are consistent with the HF 
results and support a mechanism of the soft gap working 
beyond the mean-field level. 

On the other hand, S. Chiesa et al. found only a dip 
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Fig. 10. (Color online) DOS by HF in ID at t = 1, U = 10/3, 
W = 20/3 (N B = 14) in the linear (a) and the logarithmic (b) 
plots. Electron filling is half-filled and hole-doped for fj, = U/2 
and fi = U/2 — 10/3, respectively. The straight lines in (b) are 
the scaling plots. 

of the DOS in 2D with exact diagonalization. 28 Namely, 
they found that the DOS is almost system-size indepen- 
dent within the energy resolution of the order of O.lt, 
which appears to disagree with the present result. How- 
ever, the pseudogap observed in their study may well 
correspond to high-energy part of a soft gap in our anal- 
yses. Indeed, the soft gap is restricted to very low en- 
ergies e.g., \E — Ep\ < 0.3t in our ID study. At higher 
energies, the DOS is almost system-size independent and 
looks like a pseudogap if the energy resolution becomes 
poor (e.g., the energy resolution AE = 0.2t) as shown 
in Fig. 13 being consistent with their observation. Thus 
we believe that further analyses at lower energies with 
higher energy resolution will reveal soft gaps also in 2D. 

4. Theory of Soft Hubbard Gap 

Our numerical observation of an unconventional soft 
gap urges us to explore an unconventional mechanism for 
its origin and an unprecedented low-energy excitations 
when all of the interaction, randomness and quantum ef- 
fects are combined. In § 4.1, we propose a multivallcy 
energy landscape as the origin of the soft Hubbard gap 
and also propose the resultant scaling theory. In addi- 
tion, we compare our theory with the numerical results 
in detail to test the mechanism. Although the scaling 
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Fig. 11. (Color online) Dependence of DOS on t (panel (a)), U 
(panel (b)) in ID within the HF approximation (N B = 14). The 
straight lines are the fit by power laws: (a) U = 10/3, W = 20/3, 
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Fig. 12. (Color online) (a) DOS in ID with the exact diago- 
nalization (open boundary condition): t = 1, U = 10, W = 
10, n = (7/2, (N s = 2,3,4,5,6). We average the DOS over 
3.2 X 10 7 realizations of disorder for 7V S = 6. (b) Scaling plot 
by ^(e.Mr 1 ) = Nr P f(eN<? /a ) (a = 0.075, /3 = 0.375). Shaded 
part shows the converged scaling curve. 
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theory has been reported briefly in the previous letter, 31 
we discuss in more detail for the benefit of its further ex- 
tension: Although a continuous distribution of random 
potentials is assumed for the moment for the simplicity 
of discussion, wc extend our scaling theory to discrete 
distributions of random potentials in § 4.2. Section 4.3 is 
devoted to consideration of effects of the long-range part 
of the Coulomb interaction and extension of the scaling 
to the long-range Coulomb interaction. Finally, in § 4.4, 
we summarize the scaling laws of the DOS obtained in 
this section. 



Fig. 13. (Color online) Linear plot of the same data as in Fig. 12 
with a lower energy resolution of 0.2t. 



4-1 Origin and scaling theory 

In this section, we propose the origin of the soft Hub- 
bard gap and construct its scaling theory. For simplic- 
ity without loss of generality, we restrict ourselves to a 
single-particle excitation for the electron side, namely, 
E > Ep. We consider the case of r; nt <IC £, where rj nt 
is the range of the interaction in the model; rj nt = for 
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the Anderson-Hubbard model. Even for n n t = 0, when t 
becomes nonzero, virtual hopping of electrons generates 
intersite effective interaction Uij though the zero point 
fluctuation, which exponentially decreases with the mu- 
tual distance — r 7 | as 



Uij cx exp(-6|r, - rj\), 



(4.1) 



where b is proportional to the inverse of the localiza- 
tion length. This effect is not considered in the ES the- 
ory because it regards electrons as classical particles. 
Within the HF approximation, the interaction energy 
(the second term of the Anderson-Hubbard Hamiltonian, 
cq. (1.2)) is decoupled as 



n 



U^2(\<t> m (n;V\ 2 \Mri;l)\' 



-|^m(n;4)| 2 |^(r i; t)| 2 )}, (4.2) 



where the orbital indices n and m run over all the occu- 
pied single-particle wave functions. Thus, from eq. (3.5), 
b is given by b = 2/£ within the HF approximation. 

By assuming the self-averaging of the DOS, the DOS 
averaged over the sites is obtained as 

A(E) = (f iV(Vi)AiCE,Vi)dVi) ,(4-3) 

V -oo ' {v T } 

where the symbol {Vy} denotes a set of random poten- 
tials Vi except for V\. Note that A\(E,V\) is the local 
DOS projected on the site 1 under the condition of the 
fixed Vi at the site 1. This local DOS A X (E, Vi) implic- 
itly depends on {Vy}. Here we decompose the average 
over the random potential into the two steps; namely, 
the average over V\ as described by / Pv(Vi)dVi at a 
fixed configuration of {Vy} and the subsequent average 
with respect to {Vy}. In the following, we first examine 
Vi -dependence of the single-particle excitation energy, 
namely A\{E, Vi) and obtain the Vi-averaged local DOS. 
Then we average the Vi-avcragcd local DOS over {Vy} 
to obtain the scaling of the soft Hubbard gap. 

First, we discuss Vi-dependence of Ai(E, Vi) for fixed 
{Vy}. When Vi decreases, the ground-state occupation 
of the site 1 changes from to 1 and then from 1 to 
2 at Vici and Vi C 2, respectively. A possible and typical 
ground state |<fo) at Vi > Vi c i is illustrated in Fig. 14(a), 
where the site 1 is empty. Here the total particle num- 
ber is N c = N a and the energy E (Vi). Near Vi c i but 
for V > Vid, a single-particle excited state \4>i) with 
N e = N a + 1 and the energy i?i(Vi) is defined by the 
electron configuration fixed to be the same as \cj>o) ex- 
cept for an addition of an electron at the site 1, as is 
illustrated in Fig. 14(b). In the interacting case, because 
c Ll^o) is not necessarily an eigenstate, there may be sev- 
eral eigenstates that have nonzero matrix elements with 
c i<7 \4>o), namely, single-particle excited states. However, 
for simplicity, we discuss the lowest single-particle ex- 
cited state, |0i), which dominates at low energies. One 
might think that \<pi) becomes the ground state below 
Vid, where N c = N a + 1. In this case, however, the single- 
particle excitation gap E\ — Eq vanishes at Vi c i leading 
to the absence of the gap in the Vi-averaged DOS. In- 
deed, this is what happens in the Anderson insulator at 
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Fig. 14. (Color online) Schematic illustration of (a) the ground 
state, (b) a single-particle excited state, (c) a nearly-degenerate 
state with the ground state and (d) a multiply-excited state, (e) 
Schematic of V\ dependence of excitation energies. 



(7 = 0, where the DOS exhibits no soft gap. Thus the 
numerical evidences of the soft gaps indicate that \<j>\) as 
a single-particle excited state is excluded by the electron 
correlation. 

In the ES theory, single-particle excited states are ex- 
cluded from low energies by the ground-state stability 
condition against an electron-hole excitation. This mech- 
anism, however, cannot be attributed to the formation 
of the soft Hubbard gap, because this excitonic effect is 
negligible in the case of the short-range interaction. We 
now consider the ground-state stability against excita- 
tions of more complex form. We assume a multivalley 
energy landscape, which may be characteristic to ran- 
dom systems. Namely, we assume that there exist many 
arbitrarily-low-energy excited states whose occupations 
are the same with those of |</>o) at the site 1 but whose 
configurations are globally different on other sites. In 
Fig. 14(c), we illustrate a state \(p' ) at a local mini- 
mum E' nearly degenerate with \4>q). Here the config- 
urations of \<j/ Q ) arc relaxed from \4>q) not only at the 
occupied site n nearest to the site 1 at the distance R 
but also at farther sites (> R). Figure 14(d) shows a 
single-particle excited state from \<fr' ) with the en- 
ergy E[ and the site-1 occupancy identical with \(j>i). 
Here, the two nearly-degenerate states, |0i) and \4>[) arc 
separated by a barrier, where multi-particle relaxation is 
required to reach from one to the other. Now E\ is given 
by (Vi — Ep) + Uu + Eq, where Un is the interaction 
energy between electrons on the site 1 and those on the 
site i. Note that the interaction energy |Z7ij| is nonzero 
only when the site i satisfies R < ^--R + C because 

of the localization (Keep in mind that R is the distance 
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to the nearest neighbor electrons from the site 1). On the 
other hand, E[ is given by E[ = ( Vi - E F ) + £\ U' u + E' Q , 
where U'u is again the interaction energy between elec- 
trons on the site 1 and those on the site i in the state 
\4>' ). Because Eq ~ E' and the configurations of \cj)'i) on 
these sites are different from those of \4>i), E[ is differ- 
ent from Ei typically by the amount as much as |C/ ln |. 
Now, out of many possible \4>' )s, one can choose \(f>' ) so 
that the energy E[ is lower than E\ by the amount |C/ ln |. 
Then \4>i) is indeed a state that is obtained by a multi- 
ple excitation from \4>o). This means that (</ , 'i|cf cr |</>o) = 
and \4>'i) is orthogonal to the single-particle excitations. 

Now we assume linear dependence of the excitation 
energies, Ei(Vi) and E[(Vi) as functions of Vi. Then 
E[{Vi) and E Q (Vi) crosses at Vi = V{ cl and for Vi < V{ cl 
the ground state becomes \4>'i) as illustrated in Fig. 14 
(e). Note that the excitation energy E[ — Ei is negative 
in the vicinity of Vi = V{ cl . For Vi > V{ cl , the state \<j)'i) 
is not counted in the DOS, because this state is not a 
single-particle excitation of \4>q), but rather a multiply- 
excited state. As a result, the single-particle excitation 
energy Ei—Eq has the lower bound at v\ = V{ cl . In other 
words, the energy difference A = |£a(V"/ cl ) - E[(V{ c i)\ 
should be the lowest energy of single-particle excitations 
counted in A near Vi = V[ cl in the region Vi > Vid- 
Therefore the Vi -averaged local DOS has a gap of A in 
contrast to the non-interacting case as 




P v (Vi)Ai(E, Vi)dVi oc H S {E -E F - A), (4.4) 



where H s is the Heaviside step function. The same argu- 
ment applies around Vi = Vi C 2- 

Here we make an additional comment. One might 
think that, as in the ES theory, it could be possible to 
lower the energy of \(j)i) from Ei to E[ by relaxing local 
electronic configurations only near the site n. It, however, 
always increases the energy of the electrons on the sites 
other than the site 1, because they have already been op- 
timized in the ground state. Thus a global reconstruction 
is required to lower the energy. 

Next, in order to derive the scaling of the density of 
states, we discuss the distribution function of A with 
respect to {Vy}. From the above discussion, it is reason- 
able to assume that A depends on {Vj} only through R. 
Under this assumption, from eq. (4.1), A scales as 

A(R) = aexp(-kR), (4.5) 

where a and b (oc £~ x ) are non- universal positive con- 
stants. It is clear that a is proportional to U when U is 
small, and b diverges for t — ► 0. Hereafter we neglect log- 
arithmic corrections. The localization lengths may fluc- 
tuate between the site 1 and the site n. However, this 
fluctuation of the localization lengths does not affect the 
asymptotic behavior of A(i?), eq. (4.5), because of the 
self averaging of the localization length between the site 
1 and the site n in the limit of R — > +oo. 

On the other hand, the probability distribution of R 
follows 

P(R) = a'R d - 1 exp(-b'R d ), (4.6) 



at long distances, where a' and b' are non-universal posi- 
tive constants again. Equation (4.6) means that the prob- 
ability of formation of a large void of electrons around 
the site 1 is exponentially rare. Equations (4.5) and (4.6) 
lead to 

Q(A) = P(R(A)) g 

oc (-logA) rf - 1 A- 1 exp(^-^(-logA)^ 

* A^cxp^-^-logA)^, (4.7) 

where Q(A) are the distribution function of A. Here we 
neglect the logarithmic correction term (— log A) d_1 . Be- 
cause the DOS at the energy E is proportional to the 
probability of A < \E — E F \, we obtain the scaling of the 
DOS from Eqs (4.4) and (4.7) as 

r\E-E F \ 

A(E) oc / dAQ(A) 

J 

oc eX p(-^(-log\E-E F \) d y (4.8) 

For d = 1, because the exponential and the logarithm 
functions cancel each other, this scaling reduces to a 
power law with a non-universal exponent, which is consis- 
tent with the observed power-law scaling in ID: A{E) oc 

\E — E F \ b ^ b . Non- universal power-law distributions of 
energies are common in Griffith phases. 42 Equation (4.5) 
indicates that A, namely the gap vanishes as t or U van- 
ishes because of a — > or b — ► oo. Furthermore, the 
exponent a = b'/b is expected to decrease as t becomes 
smaller because of the reduction of £. These predictions 
are consistent with our HF results in ID as shown in 
Fig. 11 (a). It should be mentioned that the power law 
for ID in the present theory is in sharp contrast with the 
ES theory, because the ES theory predicts the logarith- 
mic scaling of the DOS as eq. (3.9). For d > 1, our scaling 
predicts that the DOS vanishes faster than any power 
law, being consistent again with our HF study in 3D. In 
Fig. 15 (b), we show a scaling plot of the DOS in 3D by 
our scaling given by eq. (4.8). Indeed, the DOS fits well 
with our scaling in 3D at low energies e.g., \E—Ep\ < 0.1, 
where the DOS does not follow a power law. 

Because our scaling has explicit dimension depen- 
dence, we further test our scaling in 2D within the HF 
approximation with the Gaussian distribution of Py . Fig- 
ure 15 (b) shows the DOS with the system size lOx 10 and 
the periodic-boundary condition ati=l,[/ = 6,T / l /r = 5. 
We averaged the DOS over as many as 4.8 x 10 4 realiza- 
tions of random potentials to obtain the high-resolution 
data. As shown in Fig. 15 (b) the DOS fits well with our 
scaling obtained from eq. (4.8) rather than a power law 
for \E-E F \ < 0.5. 

In addition to the qualitative consistency between our 
scaling and the numerical results, in Fig. 16, we show a 
further numerical evidence of our theory on the quanti- 
tative level for ID. Figure 16 (a) and (b) shows A(R) 
and P(R) calculated by the following procedure: First, 
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Fig. 15. (Color online) (a) Scaling plot of DOS in 3D at param- 
eters B and C depicted in Fig. 2: B (t = 1, (7 = 4, W = 30), 
C (t = 1, [7 = 6, W = 5). We employ Lorentz broadening 
with a broadening factor 1.25 X 10~ 3 and 6.25 X 10~ 4 for A and 
B, respectively. The broken lines denote \E — Ep\ = 10~ 2 and 
10 -1 . The DOS fits well with A(E) <x exp(-(~7 log \E - E F |) 3 ) 
shown by the fitting lines for 10~ 2 < \E - E F \ < 10" 1 . (b) 
Scaling plot of DOS in 2D: t = 1, U = 6, W = 5 with the 
Gaussian distribution of Py . The inset is a linear plot of the 
same data. The black solid curves are fit by the predicted scal- 
ing, A(E) oc exp(-(7log|£ - E F \) 2 ). If A{E) followed A{E) oc 
exp(— (7 log \E — Ep\) 2 ), the data in the panel (b) would follow 
a straight line. The broken lines denotes \E — Ep\ = 0.05,0.5. 



we obtain the ground state for each realization of ran- 
dom potentials. We construct the lowest single-particle 
excited state by adding one electron to the lowest unoc- 
cupied orbital. Next we optimize the mean fields by the 
iterative scheme starting from those of the single-particle 
excited state with N e fixed. Then A is obtained as the 
difference of these two excitation energies. We calculate 
R as the distance between the center of the lowest un- 
occupied orbital, r and those of the occupied orbitals 
nearest to r in the ground state. We define the center 
of the orbital as the site that has the maximum weight. 
Fitting by eqs. (4.5) and (4.6) gives b' = 1.06 ±0.01 and 
b = 1.34±0.01. Estimated exponent ofV/b = 0.79±0.02 
is in good agreement with a = 0.85 ± 0.07 obtained di- 
rectly from the DOS as shown in Fig. 16 (c). This is a 
numerical evidence for the quantitative validity of our 
theory. 

In the above discussion, under the assumption of a 
multivalley energy landscape, we have successfully con- 
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Fig. 16. (Color online) Numerical estimates of P(R) and A(_R) 
at t = 1, U = 10/3, W = 20/3 and fj, = U/2 - 10/3. Fittings by 
eqs. (4.5) and (4.6) as are shown in (a) and (b) by dotted lines 
give b' = 1.06 ± 0.01 and b = 1.34 ± 0.01. The ratio b' jb (dashed 
lines in (c)) well agrees with a obtained by the fitting (solid line) 
of the DOS in (c). 



structed the scaling theory which is consistent with all 
the numerical results. In the presence of a multivalley 
energy landscape, there are many excited states whose 
configurations are globally different from those of the 
ground state, which is a non-trivial combined effect of 
the electron correlation and randomness. Indeed, in the 
field of spin glasses, it is known that a multivalley en- 
ergy landscape emerges concurrently with replica sym- 
metry breaking (RSB) within a spin-glass phase. Exten- 
sive studies on replica symmetry breaking (RSB) in fi- 
nite dimensions have been carried out especially in clas- 
sical Ising models (Ising Edward- Anderson model). 43 On 
the other hand, our Hartree-Fock results indicate that a 
multivalley energy landscape exists over the whole in- 
sulating phases for d = 1, 2 and 3. Furthermore, al- 
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though the cxact-diagonalization results arc restricted to 
the strongly-localized region in ID because of the severe 
limitation of system sizes, the exact-diagonalization re- 
sults indicate robustness of the multivalley energy land- 
scape against quantum fluctuations beyond the mean- 
field level. Further theoretical studies along this line are 
intriguing future subjects. 

Our conclusion disagrees with the DMFT results 23 as 
well as several other mean-field studies 25,30 which sup- 
port absence of the soft gaps. This may be because they 
ignore spatial correlations. The latter ignores inhomo- 
geneity of the electronic structures. Indeed, a DMFT 
study improved by partially taking account of the in- 
tersite self-energy retrieves the suppression of the DOS 
near Ef to some extent. 29 

In contrast to our scaling, the power law was proposed 
to interpret the photoemission experiments. 15, 16 How- 
ever, as shown in Figs. 9 and 15, the asymptotic be- 
havior of the DOS is restricted to low energies, namely 
\E — E-p\ < O.lt and high-energy part of the DOS seems 
as if it approximately follows a power law with a non- 
universal exponent in our 3D HF study. Because hop- 
ping integrals between d orbitals on nearest-neighbor Ru 
atoms are on the order of 0.1 cV for SrRu03, 44 the 
asymptotic behavior of the soft Hubbard gap may be 
restricted to the energy region lower than 10 meV. Thus 
photoemission experiments with the energy resolution on 
the order of 1 meV are desired to observe the present 
asymptotic behavior clearly. Although this level of reso- 
lution has become possible recently, such high resolution 
has not been utilized so far for the present purpose in the 
literature. 15, 16 We believe that our paper provides incen- 
tive for such high-resolution photoemission experiments. 
In addition to the high-resolution photoemission, other 
experiments accessible to low energies arc highly desired 
for experimental confirmation of the present theory. For 
example, the DC transport measurement is suitable for 
investigating the density of states in insulating phases at 
low energies i.e. T < 300 K {k B T < 30 meV). Actually, 
in § 5, we discuss the temperature dependence of the DC 
resistivity in the presence of the soft Hubbard gap. 

4-2 Extension of scaling theory to discrete distributions 
of random potentials 
Although we assume a continuous distribution Py in 
the previous sections, the obtained scaling, cq. (4.8) does 
not depend on the distribution. Even for a discrete dis- 
tribution of a binary alloy form, one obtains the same 
scaling as that for continuous distributions by modifying 
the above discussion slightly. In the above discussion, 
we first average the single-particle excitation spectrum 
over Vi at fixed configurations of {Vj}, which is a set 
of random potentials Vi except for V\ . For a continuous 
distribution, as described in eq. (4.4), the single-particle 
excitation energy distributes continuously above a lower 
bound A(i?) depending on {Vj} through R. Here R is 
the minimum distance from the site 1 to the occupied 
sites in the ground state. This continuous distribution is 
a key point for obtaining our scaling. For a discrete dis- 
tribution, however, the single-particle excitation energy 
distributes discretely above A(i?), because V\ is discrete. 
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Thus we need an extension of our consideration here. We 
further average the single-particle excitation spectrum 
over {Vj} with a certain common value of R. Then the 
distribution of the single-particle excitation energy be- 
comes continuous above A(i?), because of fluctuations of 
the electronic structures at the distance R from the site 
1 in the ground states. After a further average over R, 
one obtains the same scaling as that for the continuous 
distribution. 

We show a numerical evidence in ID within the 
Hartrce-Fock approximation in Fig. 17. We employ the 
bimodal distribution: Vi = ±W with the symmetric 
probability P(W) = 1/2 and P{-W) = 1/2. Although 
the density of states has a complex peak structure due to 
the discreteness of the bimodal distribution, the density 
of states clearly exhibits a soft gap near the Fermi energy 
as shown in Fig. 17 (a). Indeed, as shown in Fig. 17 (b), 
the soft gap follows a power law, which is consistent with 
the above discussion. 

(a) 0.7 fl ' ' ' 11 
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Fig. 17. (Color online) Density of states in ID within the Hartrec- 
Fock approximation with 24 sites: t = 1,17 = 2,W = 2. (a) 
Linear plot, (b) Double logarithmic plot. The solid line in (b) is 
the fit with a power law. 



4-3 Effects of the long-range part of the Coulomb inter- 
action 

Even near metal-insulator transitions where the 
screening is strong, the Coulomb interaction remains long 
ranged in insulating phases, though its amplitude may be 
small. Thus the scaling of the DOS deviates from eq. (4.8) 
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in real materials, as the long-range part becomes domi- 
nant at low energies. Note that eq. (4.8) is the asymp- 
totic scaling for the model with short-range interaction 
only. In this section, we discuss effects of the remaining 
long-range part of the Coulomb interaction. 

Owing to the long-range part of the Coulomb inter- 
action, Efros et al. showed that a ground-state stability 
against an exciton excludes low-energy single-particle ex- 
citations and generates power-law soft gaps as A{E) oc 
^-E^- 1 (sec cq. (l.l)). 6 However, since they consid- 
ered the stability condition only against electron-hole ex- 
citations, the DOS may vanish faster than the ES scaling 
at lower energies because of stability conditions against 
further many-particle excitations. Indeed, Efros obtained 
a correction of the scaling in 3D by considering multi- 
clcctron-hole excitations as 45 



A{E) cx exp(-(eoA) 1/2 ), 



(4.9) 



where e = \E — Ep\ and eo is a constant. For d — 1, 2, he 
claimed that there is no further correction in ES scaling. 

However, since our scaling even for the short-range in- 
teraction shows faster decay of the DOS than any power 
law in the presence of a multivallcy structure, it is nat- 
ural to infer that the DOS should decay even faster in 
the presence of the long-range interaction and certainly 
beyond the ES theory at low energies instead of eq. (4.9). 
Thus we extend our scaling theory to the case of the long- 
range Coulomb interaction. For the long-range Coulomb 
interaction, eqs. (4.5) and (4.6) are modified as 

A(R) = aR'\ (4.10) 

a'i? d - 1 cxp(-6'i? d ). 



P{R) 

From eqs. (4.10) and (4.11), we obtain 
Q(A) = P(R) 1 



(4.11) 



dA 



* d a'A- d - 1 exp(-a < VA- d ). (4.12) 

In a manner similar to the short-range interaction 
(eq. (4.8)), the scaling of the DOS is given by 

A(E) oc / Q(A)dA 



x 



/ A- d - 1 exp(-a d 6'A- d )dA 
Jo 



1 



-exp(-a d 6'e- d ) 



a d b'd 
oc exp(-a d b'e~ d ). 



(4.13) 



Indeed, this scaling shows faster decay of the DOS than 
eq. (4.9) proposed by Efros for 3D. Moreover, eq. (4.13) 
reveals that the ES scaling must be modified even for 
d=l,2. 

4-4 Summary of scaling laws 

In Table I, we summarize the scaling laws of the DOS 
in four kinds of models: (A)/(B) short-range interaction 
without/with a multivalley energy landscape, (C)/(D) 
long-range interaction without/ with a multivalley energy 
landscape. Now we consider how the scaling of the soft 



gap depends on the energy for the case (D), which cor- 
responds to realistic materials. At energies higher than 
the energy scale of the long-range part, the formation 
process of the soft gap is dominated by the short-range 
part. Thus, in this energy region, the soft gap follows the 
scaling of the soft Hubbard gap (eq. (4.8)). As the long- 
range part of the Coulomb interaction becomes dominat- 
ing at low energies, the scaling of the DOS crosses over 
from eq. (4.8) to the ES scaling eq. (1.1) and further to 
eq. (4.13). It should be mentioned that the intermediate 
scaling laws may not always be observed clearly. 

5. Transport Properties 

In this chapter, we discuss the temperature depen- 
dence of the DC resistivity in the presence of the soft 
Hubbard gap or the soft Coulomb gap. Although the 
transport coefficient is determined by the two-particle 
(electron-hole) correlations and is not necessarily identi- 
cal with the single-particle excitations measured in the 
DOS, the DC transport measurement is a useful and 
good tool for investigating the DOS at low energies with 
high resolution in insulating phases if the transport is 
dominated by independent single-particle excitations. 

5.1 Temperature dependence of the DC resistivity 
5.1.1 Insulators with a hard gap 

When a gap is still open for weak disorder, thermally- 
excited carriers dominate the conduction and the tem- 
perature dependence of the DC resistivity follows the 
Arrhcnius law as 

'TV 
, T 

where To is a constant. 



Po cxp 



(5.1) 



5.1.2 Anderson insulators without soft gaps 

If a gap is closed for stronger disorder, localized states 
near Ep dominates the conduction. Mott showed that 
at sufficiently low temperatures, conduction results from 
electron hopping between localized state within a nar- 
row band near Ef, which is called variable-range hop- 
ping (VRH). 46 He also showed that provided there is a 
non- vanishing DOS at Ef, the temperature dependence 
of the DC resistivity exhibits universal behavior as 

i/(d+i) 



P = Po exp 



(5.2) 



5.1.3 Anderson insulators with soft gaps 

Because VRH explicitly depends on the DOS near Ef , 
soft gaps modify the temperature dependence of the re- 
sistivity qualitatively. Now we discuss scaling of the DC 
resistivity in the presence of a soft gap. As discussed in 
§ 4.3, the scaling of the DOS crosses over from eq. (4.8) 
to the ES scaling eq. (1.1) and further to eq. (4.13) as 
the long-range part of the Coulomb interaction becomes 
dominating at low energies. With this evolution of the 
crossover, p also crossovers as a function of the temper- 
ature. In the following, we discuss scaling of the DC re- 
sistivity for each energy region separately. 

(1) Energy regions (D.2) and (D.3) in Table II 
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(D.l) exp ( 


-7(-loge) d ) 


eq. (4.8) 


(HEs) 


(D) Long-range Coulomb interaction 


(D.2) : 




eq. (1.1) 


i 


with multivalley energy landscape 


(D.3) exp(- 


-(eoA) 1/2 ) (3D) 


eq. (4.9) 


1 




(D.4) exp(- 


-/3e~ d ) 


eq. (4.13) 


(LEs) 



Tabic I. Summary of scaling laws of DOS four kinds of models: (A)/(B) short-range interaction without/with a multivalley energy 
landscape, (C)/(D) long-range interaction without/with a multivalley energy landscape. Shaded part denotes the novel scaling laws 
obtained in this paper. The scaling of the case (C.2) is expected to appear in energies lower than that of the case (C.l). The case (D.l) 
corresponds to the energy scale where the short-range part of the Coulomb interaction is dominant, while the cases (D.2-4) correspond 
to that where the long-range part is dominant. In the cases of the long-range interaction, scaling crossovers arc expected as a function 
of the energy. Abbreviations are: High energies (HEs); Low energies (LEs). 



Efros et al. 6 showed that in the presence of the soft 
Coulomb gap with the power law A{E) oc \E — _E F | d_1 
(eq. (1.1)), the temperature dependence of the DC resis- 
tivity is modified regardless of the spatial dimension d 
as 

, 1/2] 

p = po exp 1 



(5.3) 



Even in the case (D.3) of Table II where the DOS is scaled 
as exp(— (eo/e) 1//2 ), Efros concluded that the DC resis- 
tivity follows eq. (5.3), because the excitation spectrum 
of the particle screened by excitons that determines the 
DC transport still follows the power law \E — Ep^ 1 . 45 

(2) Energy region (D.l) in Table II 

On the other hand, at high energies where the short- 
range part of the Coulomb interaction is dominant, the 
DOS follows the scaling of the soft Hubbard gap eq. (4.8). 
Here we discuss the DC transport in the presence of 
the soft Hubbard gap without considering the long-range 
Coulomb interaction. First, we assume the DOS in the 
form 



A( e )=a|e|^exp(-( 7 |loge|) rf ), 



(5.4) 



where e = \E— Ep\. The power-law correction term |e|' 3 is 
a slight generalization from eq. (4.8). Under this assump- 
tion, we obtain the resistivity up to the leading order at 
low temperatures via variable-range hopping as 



po exp c 



exp(-c 1 |log(fc B T)| 1 /rf) 



k B T 
-l/d 



(5.5) 



for d > 1, Co = l+a-^-^fj)- 1 ^, a = £/2 and a > 0. For 
details of the derivation of eq. (5.5), readers are referred 
to Appendix A.l. 

(3) Energy region (D.4) in Table II 

We next discuss the modification of the DC resistivity 
at energies lower than those justified by the ES scal- 
ing, by starting from eq. (4.13). We derive the temper- 
ature dependence of the DC resistivity in a way similar 
to the case of the soft Hubbard gap. Following a pro- 
cedure similar to the case of the soft Hubbard gap (see 



Appendix A. 2), we obtain the resistivity up to the lead- 
ing order as 

p = pocxpl — I. (5.6) 

The scaling laws of the DC resistivity obtained in the 
above discussion are summarized in Tabic II. In Fig. 18 
(a), we compare the obtained DC resistivity for the mod- 
ified ES scaling, namely eq. (5.6) with that for the ES 
scaling, namely eq. (5.3). Indeed, the DC resistivity for 
the modified ES scaling diverges at low temperatures 
faster than that for the ES scaling. However, as shown in 
Fig. 18 (b), the DC resistivity for the modified ES scaling 
is almost indistinguishable from the Arrhenius law just 
by tuning the activation energy. 

Here we propose a procedure to distinguish the modi- 
fied ES scaling from the Arrhenius law clearly: 

(1) One estimates po in the scaling for the modified ES 
scaling (eq. (5.6)): From eq. (5.6), one obtains 

{P/df/ d \\og(k B T)\- l / d 



log(p) = log(po) 



log(po) (1/T-. +0). 



(5.7) 



Thus one can estimate po easily by a linear extrap- 
olation of log(p) to the limit of 1/T — > +0 with 
respect to 1/T. 

(2) Next, one plots the logarithm of normalized re- 
sistivity multiplied by the temperature T\og{p/po) 
against | log(T)| _1 / d as shown in Fig. 18 (c). In this 
plot, log(p/po) is proportional to |log(T)| _1 / d for 
the modified ES scaling, while the normalized resis- 
tivity \og(p/po) is constant when the experimental 
data follows the Arrhenius law. These two functions 
can be distinguished clearly from each other. 

We believe that the validity of the present theory can be 
tested against experiments by using this plot provided 
that po does not contain other extrinsic factors of tem- 
perature dependence out of the present consideration, 
such as the volume expansion. 
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Fig. 18. (Color online) Comparison of the DC resistivity for the 
modified ES scaling, eq. (5.6) with (a) that for the ES scaling, 
and (b) the Arrhenius law (the dashed and dotted lines). We 
note that the modified ES scaling shown by the solid line is 
hardly distinguishable from the Arrhenius law (the dotted line) 
by fitting the activation energy. In eq. (5.6), we take fee = 1 
and (ft/d) 1 ^ = 1. (c) In this plot, all the data which follow 
the Arrhenius law converge to a constant, namely unity irre- 
spective of the activation energy, while log(p/po) is proportional 
to | log(T)| _1 / ci for the modified ES scaling. This distinguishes 
these two scaling laws clearly. 



5.2 Comparison with experiments 

In this section, wc compare our formulae for the DC 
resistivity, eq. (5.5) and eq. (5.6) with transport proper- 
ties of SrRui-^TiajOs, 14 where the breakdown of the ES 
scaling was indicated in photoemission experiments. 15 ' 16 
In Fig. 19, we show the experimental data of the DC re- 
sistivity of SrRui^Ti^Os at x = 0.6, where the ground 
state is insulating. The experimental data fit well with 
the ES scaling for k B T < 0.005 eV (T < 60 K), indicating 
that the ground state is an Anderson insulator. At higher 
temperatures (0.005 eV < k B T < 0.011 eV, 60 K < T < 
130 K), however, they fit well with the scaling obtained in 
the presence of the soft Hubbard gap for 3D (eq. (A-2)), 
being consistent with our theory. The fitting parameters 
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Fig. 19. (Color online) Scaling analyses of the DC resistivity of 
SrRui—jjTizOa (the experimental data illustrated as filled circles 
are from Ref. 14). (a) The logarithmic plot. The panels (b), (c) 
and (d) show scaling plots for the ES scaling (T < 55K), the 
soft Hubbard gap (55K < T < 135K) and the Arrhenius law 
(135K < T), respectively. The vertical dotted lines denote T = 
55 K, 135 K. 
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are as follows: log(p ) = 9.0 ± 0.2, log(co) = 0.94 ± 0.7 
and ci = 2.6 ±0.4. On the other hand, for higher temper- 
atures (0.011 eV < k B T), the experimental data deviate 
from the scaling for the soft Hubbard gap, being con- 
sistent with the our 3D HF study where the asymptotic 
behavior of the soft Hubbard gap is restricted to low en- 
ergies (< O.li) as shown in Figs. 9 and 15. On the other 
hand, we cannot find the region where the modified ES 
scaling predicted by eq. (5.6) is satisfied even in the low- 
est temperature available in the experiment (~ 30 K). 
We believe that further experiments at lower tempera- 
ture will reveal the existence of the correction to the ES 
scaling. 

Another interpretation of the experimental data for 
0.005 eV < k B T < 0.011 eV may be that the ES scal- 
ing (k B T < 0.005 eV) and the Arrhenius law (k B T > 
0.011 eV) cross over there; the Arrhenius law arises from 
the hopping conduction between nearest-neighbor local- 
ized states, which dominates the DC resistivity at high 
temperatures instead of the variable-range hopping. In- 
deed, for 0.011 eV < k B T < 0.026 eV (130 K < T < 
300 K), they seem to follow the Arrhenius law with a 
small activation energy 0.030 ±0.001 eV. In order to ex- 
clude or support the possibility of this crossover, further 
experiments as functions of composition x in the vicinity 
of the metal-insulator transition are desired. 

6. Summary and discussion 

In summary, we have examined the ground state and 
single-particle excitations of the 3D Anderson-Hubbard 
model within the Hartree-Fock approximation. In § 3, 
although only the short-range interaction is present in 
the model in contrast to the Efros-Shklovskii theory, we 
have found an unconventional soft gap over the whole 
insulating phases. Namely, we have found that the DOS 
vanishes toward the Fermi energy faster than any power 
law indicating the formation of an unconventional soft 
gaps, which we call soft Hubbard gap. Further numerical 
evidences to support this unconventional soft gap has 
been given within the Hartree-Fock approximation, and 
further with the exact diagonalization in ID. In contrast 
to the 3D case, a power-law scaling is satisfied in ID. 

In § 4, based on the picture of a multivalley energy 
landscape, we have constructed a phcnomenological the- 
ory, being qualitatively as well as quantitatively consis- 
tent with the numerically-obtained scaling of the soft 
Hubbard gap. There, further support for the scaling the- 
ory has been given in 2D within the Hartree-Fock ap- 
proximation. Moreover, by considering effects of the long- 
range part of the Coulomb interaction remaining un- 
screened at low energies, we have obtained a novel scaling 
of the soft gap beyond the Efros-Shklovskii theory. This 
novel scaling is caused by a multivalley energy landscape. 
Finally, we have proposed scaling crossovers of the DOS 
as the long-range part becomes dominant at low tem- 
peratures. Scaling laws and their expected crossovers are 
summarized in Table II. 

In § 5, we have derived the temperature dependence of 
the DC resistivity in the presence of the soft gaps. Pos- 
sible experiments to verify the present theory has been 
proposed. The scaling laws of the density of states and 



the DC resistivity are summarized in Table II. Finally, we 
have shown that the present theory is indeed consistent 
with the experimental data for SrRui_ x Ti x 03. 

In order to readdress the essence of the present theory 
and stimulate further theoretical studies, we now discuss 
possible extension of the present theory. In § 4, we have 
shown that a gapless collective mode accompanying a 
multivalley energy landscape causes relaxations from a 
single-particle excited states to a multiply-excited states 
not counted in the single-particle DOS, which directly 
leads to the formation of the unconventional soft gap. In 
other words, the single-particle excited states are elim- 
inated from low energies by the ground-state stability 
against the collective excitations at low energies. How- 
ever, there are actually many other kinds of gapless col- 
lective modes. A typical example is a spin wave accom- 
panying magnetic ordering with spontaneous symmetry 
breaking, which is not included in the present Hartree- 
Fock approximation. If these collective modes also cause 
the relaxation from a single-particle excited states to a 
multiply-excited states, the present scaling of the den- 
sity of states will equally be valid, which extends the 
applicable scope of the present theory Whether or not 
the present theory can be extended to general collective 
modes is left for a future challenge. 

Recently magnetic hard gaps have been observed in a 
large number of disordered materials such as In-doped 
CdMnTe, 48 - 49 amorphous GeCr-films, 50 doped Si 51 and 
amorphous Sii-^Mn^-films 52 by measuring the temper- 
ature dependence of the DC resistivity. These gaps were 
identified to be of magnetic origin by reduction of the gap 
sizes by an increased external magnetic field. Because the 
DOS vanishes toward the Fermi energy even faster than 
a power law in the presence of the soft Hubbard gap, the 
soft Hubbard gap causes the rapid divergence of the DC 
resistivity slightly slower than the Arrhenius law at low 
temperatures as summarized in Table II. Comparisons of 
the present theory with these experiments are interesting 
issues. 
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Appendix: Derivation of temperature depen- 
dence of DC resistivity 

In this appendix, we derive the temperature depen- 
dence of the DC resistivity in the presence of the 
soft Hubbard gap or the soft Coulomb gap within the 
variable-range hopping. 

A.l Soft Hubbard gap: Energy region (D.l) in Table II 
First, we derive the DC transport in the presence of 
the soft Hubbard gap. We assume the DOS in the form 

A( J B)=a| e |^exp(-( 7 |loge|) d ), (A-l) 
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Model 



Scaling of DOS 



DC resistivity 



Temperature 



Clean Mott insulator 
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exp(T /r) 

exp((T /T)W+U) 



(A) Short-range interaction 

without multivalley energy landscape 



A >0 



r-n (r . exp(- ci | logCfce^l^T T 

exp ^c kBT j 



(B) Short-range interaction 

with multivalley energy landscape 



exp (-7(-loge) d ) 
eq. (4.8) 



eq. (5.5) 



e dLl exp ((T /T)V*) 

exp(-( £o /6) 1 / 2 ) (3D) exp((r /r) 1 / 2 ) 

exp(-7(-lo K e) d ) exp fc n Gxp( " Cl "; 
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(D) Long-range Coulomb interaction 
with multivalley energy landscape 
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• • eq. (5.5) 
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(D.4) exp(-/3e" d ) 


exp 
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■•■ eq. (4.13) 




■ ■ eq. (5.6) 
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Table II. Summary of scaling laws of DOS and DC resistivity for four kinds of models: (A)/(B) short-range interaction without/with 
a multivalley energy landscape, (C)/(D) long-range interaction without/with a multivalley energy landscape. Shaded part denotes the 
novel scaling laws obtained in this paper. The scaling of the case (C.2) is expected to appear in energies lower than that of the case 
(C.l). The case (D.l) corresponds to the energy scale where the short-range part of the Coulomb interaction is dominant, while the 
cases (D.2-4) correspond to that where the long-range part is dominant. In the cases of the long-range interaction, scaling crossovers 
are expected as a function of the temperature. Abbreviations are: High temperatures (HTs); Low temperatures (LTs). 



where e = \E—Ep\. The power-law correction term |e|' 3 is 
a slight generalization from eq. (4.8). In the following, we 
discuss p for the soft Hubbard gap without considering 
the long-range Coulomb interaction. 

In order to determine the energy window £q dominat- 
ing the transport for a given temperature T, we minimize 
the resistivity 



p = p Q exp 



1 



(A-2) 



N 1 / d (e )a k B T J ' 

derived from the hopping transport within the energy 
window of the width eo, where a = £/2. We define the 
number of electrons within the energy window, N(eq) as 

A^(eo) = / A(E)dE 



-€Q-\-E-p 

" + 00 

exp(- 7 V - (p + l)x)dx 



2a 



2a 



P+l 



T exp(- 7 V-(/?+l) C ) 
e^ +1 cxp(- 7 d |logeo| d ), (A-3) 



which is justified for c 3> 1. In the fourth line, we expand 
the incomplete gamma function Y(z,p) in p as 

T(z,p) 



jl + £^-l)(^-2)---(z-n)J(A.5) 



where x = — log(e), c = — log(eo)- Here wc employ an 
approximation 

/ exp(-ax d - x)dx = / ^—{-a) n {x d ) n e- x dx 
Jc Jc n=0 n. 

00 

n! 

n=0 

= exp(-ac d -c), (A-4) 



By differentiating eq. (A-2) with respect to e , we ob- 
tain 

dp 
de 

Thus the minimization condition of the resistivity is 
given by 

2 ^°) = J_ (A .7) 
adiVVd+i(e ) fc B T' K 1 

By taking logarithms of both sides of eq. (A- 7) and sub- 
stituting eqs. (A-l) and (A-3), we obtain 

7 d |lo geo | d + (/3 + d+l)|logeo| =d|Iog(feT)|, (A-8) 

up to the next leading term. In general, a solution of an 
equation for x, 

a a x + b x 1/d = c (A-9) 

is given by 

-1 , l/d -1-1/d . ^.y -1+2/ds t a 1n \ 

x = c a -b c ' a + O{c ). (A-10) 

Thus, for d > 1, we obtain the solution of eq. (A-8) up 
to the next leading term in the low-temperature limit as 

7 d |loge |' i = dX-{l3 + d+l)d l ' d 1 ~ l X 1,d 

+0{X- 1+2 / d ), (A-ll) 



expf-dV^-ixVd 

+(p + d + i)d- 2+2 / d 1 - 2 x- l+2 ' d 

+0(X- 2+3 / d )), (A-12) 
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where X represents |log(fcBT)|. By substituting these 
equations into eq. (A-3), we obtain 

iV(eo) = ■^ T (k B T) d cMd 1/d+1 l~ 1 X 1 / d ). (A-13) 

Thus the first term in the exponential function in 
eq. (A-2) up to the leading term is given by 

1 exp(-d 1 / d ~/- 1 X 1 / d ) 



(A14) 



On the other hand, the second term in the exponential 
function up to the leading term in eq. (A-2) is given by 

e cxp(~d 1 / d 7 - 1 A- 1 / d ) 
k B l k B l 

Thus we obtain the resistivity up to the leading order at 
low temperatures as 

expt-dllog^BT)! 1 ^)' 



P = Po exp c 



k B T 



(A-16) 



for d > 1, c = 1 + a _1 (^r) _1/d and ci > 0. Equa- 
tion (5.5) gives the scaling of the resistivity for our soft 
Hubbard gap (given by eq. (4.8)) obtained for the model 
with the short-range interaction only. 

A. 2 Modified ES scaling: Energy region (D.4-) in Ta- 
ble II 

We next discuss the modification of the DC resistivity 
at energies lower than those justified by the ES scaling, 
by starting from eq. (4.13). We derive the temperature 
dependence of the DC resistivity in a way similar to the 
case of the soft Hubbard gap. 

We integrate the DOS with respect to energy using the 
expansion of the incomplete gamma function and obtain 

N(e ) = / A(E)dE 



-c+E F 



1 



= 2al3- 1 d~ 1 T(--,^ d ) 

~ 2a l 3- 1 d- 1 e d+1 exp(-/3|e |- d ). (A- 17) 

By taking logarithms of both sides of eq. (A-7), we 
obtain 

/3cT V + (1 + 1/d) 2 log( eo " d ) = | log(fc B T)|. (A18) 

up to the next leading term. In general, the solution of 
an equation for x in the limit of Co — > +oo, 

a x + b log(x) = c (A- 19) 

is given by 

x = c a ' 1 -a ( 7 1 feolog(co) + 0(c°). (A-20) 

Thus we obtain the solution of eq. (A- 18) up to the next 
leading term in the low-temperature limit as 

e d = dp~ x X -dp" 1 (1 + 1/d) 2 log(X) 

+O(X ), (A-21) 

e = {p/d) 1 / d X- 1 l d 

+ (P/d) 1/d d- 1 (l + l/d) 2 X-( 1+1 '^ log(X) 

+0(X- {1+1/d) ), (A-22) 



where X represents | log{k B T)\ again. Thus the first term 
in the exponential function in eq. (A-2) up to the leading 
term is given by 

i |iog(fc B r)|-(W*) 

{C/a) — ,(A-23) 



N(e y/ d a k B T 

where C is a constant. On the other hand, the second 
term in the exponential function up to the leading term 
in eq. (A-2) is given by 

jo_ (/W /d |iog(fc B r)r 1/d , . 

k B T k B T K 1 

Thus we obtain the resistivity up to the leading order as 
f(P/d) 1 / d \\og(k B T)\- 1 / d ' 



Po cxp 



k B T 



.(A-25) 
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